
# This cell will import packages that we will need to run the scopes and collect and analyze data.

import pyvisa as visa  # This is python package that controls the usb communications
import numpy as np # This package will allow you to build arrays of data points 
import time # This is sometimes necessary to pause the code to allow a command to complete on one of the instruments.
import matplotlib.pyplot as plt # This allows you to plot data on graphs
from scale_scope import *
from scipy.fft import fft, ifft, fftfreq
import do_spec_2
#import Transfer_Function_Code as TFC

from os.path import exists # use this to protect from overwriting an existing data file.
import datetime # use to create unique data file names


"""
This cell will set up the measurement file for future analysis. This file will include a full accounting of the 
experimental setup. This will include:
    Date of measurement in "mo/da/year" format
    Experimentor "Your Name Here"
    Oscilloscope used (usually RIGOL DS1054Z)
    Function Generator used (usually RIGOL DG1022Z)
    Audio amplifier on the output signal (Yes or No)
    Audio amplifier Channel A settings ( Volume, Presence)
    Porous Sample used (consult list of disks in separate file)
    Air backing thickness (in cm)
    Microphone A location (in cm, relative to sample surface)
    Microphone B location (in cm, relative to sample surface)
    Number of frequencies in data run
    Minimum Frequency
    Maximum Frequency
"""
# this is the full path to the data directory, it is not stored in the Github directory as there will be large
# amounts going forward which should not be uploaded to Github
fpath = "/home/research/Desktop/NewScope/Tube_Data/"
# filename for the data file to be recorded. it should folow this format: mo_da_year_starttime.dat
td = datetime.date.today()
ts = datetime.datetime.now()

#"whotheplugbelongsto-Areaofporosity-numberofpores-airbackingringheight(cm)-microphoneorder
fname = "DrFreds-A50-900-2-AB_test" #AB=channel #1 microphone is on top and channel 2 is on bottom
path_to_file_dat = fpath + fname + ".dat"
file_exists = exists(path_to_file_dat)
if file_exists:
    print("old file exists, rename output file")
else:
    print(path_to_file_dat)
fdat = open(path_to_file_dat, "w")
print("Data file is: "+path_to_file_dat)
path_to_file = fpath + fname + ".info"
finfo = open(path_to_file, "w")
print("Info file is: "+path_to_file)
#path_to_file = fpath + fname + ".freq"
#ffreq = open(path_to_file,"w")
print("file open for writing")
finfo.write(str(td.year)+"-"+str(td.month)+"-"+str(td.day)+"\n") #Date of measurement
finfo.write("Cole\n") # responsible person
finfo.write("RIGOL DS1054Z\n") # Scope used
finfo.write("RIGOL DG1022Z\n") # Function Generator used
finfo.write("Yes\n") # Amplifier used yes or no
finfo.write("1, 0\n") # Amplifier settings channel A (if not used pu NA)
finfo.write("Sample A50-900\n") # porous sample used (see list of disks)
finfo.write("2.0 cm\n") # Air backing thickness (in cm)
finfo.write("11.0 cm\n") # Microphone A location (in cm, relative to sample surface)
finfo.write("2.75 cm\n") # Microphone B location (in cm, relative to sample surface)
df = 10.0 # Measurement frequency step size in Hz
finfo.write(str(df) + " Hz\n") # frequency step used in data run
fmin = 100 # Hz
finfo.write(str(fmin) + " Hz\n") # Minimum Frequency
fmax = 2000.0 # Hz
finfo.write(str(fmax) + " Hz\n") # Maximum Frequency
finfo.close()

/home/research/Desktop/NewScope/Tube_Data/DrFreds-A50-900-2-AB_test.dat
Data file is: /home/research/Desktop/NewScope/Tube_Data/DrFreds-A50-900-2-AB_test.dat
Info file is: /home/research/Desktop/NewScope/Tube_Data/DrFreds-A50-900-2-AB_test.info
file open for writing


# This cell will initialize the control of the usb ports and indentify any scopes or function generators that
# are attached.

rm = visa.ResourceManager('@py') #use pyvisa as the backend and define a resource manager
instruments = rm.list_resources() # get a list of instrument attached to the usb ports

# print out all of the instruments attached
print('A list of all instruments attached to the USB ports')
print(instruments)
print(' ')

# filter this list for Rigol scopes
Rigol_scopes = list(filter(lambda x: 'DS1' in x or 'HDO1' in x, instruments)) 

# filter this list for Rigol funtion generators
Rigol_funcgen = list(filter(lambda x: 'DG1' in x, instruments)) 

# print out the list of Rigol scopes
print('A list of the Rigol scopes attached to the USB ports')
print(Rigol_scopes)
print(np.shape(Rigol_scopes)[0])

# print out the list of Rigol funtion generators
print('A list of the Rigol function generators attached to the USB ports')
print(Rigol_funcgen)

A list of all instruments attached to the USB ports
('USB0::6833::1602::DG1ZA224203460::0::INSTR', 'USB0::6833::1230::DS1ZA224211060::0::INSTR')
 
A list of the Rigol scopes attached to the USB ports
['USB0::6833::1230::DS1ZA224211060::0::INSTR']
1
A list of the Rigol function generators attached to the USB ports
['USB0::6833::1602::DG1ZA224203460::0::INSTR']


# This cell will initialize the scopes and function generators so that we can pass commands and receive data back 
# from them. Each instrument is initially set to its factory settings as a precaution. We will later put in the 
# specific settings that we want to use for our measurements. By starting from the factory settings, we can be 
# sure that there are no settings left over from a previous user.

print("The number of scopes is: ", np.shape(Rigol_scopes)[0])
print("The number of fucntion generators is: ", np.shape(Rigol_funcgen)[0])

myscopes = [] # allow from more than 1 scope to be controlled
# Initialize the scope
for i in range(np.shape(Rigol_scopes)[0]):
    myscopes.append(rm.open_resource(Rigol_scopes[i], timeout=100000, chunk_size=1024000))
    myscopes[i].write(":*RST") # This is the command to set the instrument to factory default settings
    print("scope ", i, " found")

# select the scope to use for data acquisition
scope = myscopes[0]

# Initialize the fucntion generator
for i in range(np.shape(Rigol_funcgen)[0]):
    funcgen = rm.open_resource(Rigol_funcgen[i], timeout=100000, chunk_size=1024000)
    funcgen.write(":*RST") # This is the command to set the instrument to factory default settings
    print("function generator found")

# give each instrument some time to reset before you go in and start changing parameters
time.sleep(5) # Wait 5 seconds to continue

The number of scopes is:  1
The number of fucntion generators is:  1
scope  0  found
function generator found


# This cell will set the scope up to make the measurements that we are interested in. It will set the channel
# we are looking at as well as the timebase and scale of the channel.

# We need to use both channels 1 and 2. 

# Channel 1 will connect to mic 1:
scope.write(":CHAN1:DISP ON") # Turn off the display of channel 1
scope.write(":CHAN1:PROB 1") # Set the probe to a x1 probe
scope.write(":CHAN1:COUP AC") # set the couplong on channel 1 to DC
scope.write(":CHAN1:SCAL 0.2") #Set the scale to 0.2 V per division (square)

# Channel 2 will be the trigger channel:
scope.write(":CHAN2:DISP ON") # Turn off the display of channel 1
scope.write(":CHAN2:PROB 1") # Set the probe to a x1 probe
scope.write(":CHAN2:COUP AC") # Set the coupling on channel 2 to DC
scope.write(":CHAN2:SCAL 2.0") # Set the scale to 2.0 V per division (square)

"""
    Need to revamp the process for collecting data as there are only 2 channels. Should be able
    to bring a trigger signal into the back to control phase for the recorded signals. Currently
    the only TRIG command that is valid is the sweep command. 
"""


# Channel 3 will connect to mic 2:
scope.write(":CHAN3:DISP ON") # Turn off the display of channel 1
scope.write(":CHAN3:PROB 1") # Set the probe to a x1 probe
scope.write(":CHAN3:COUP AC") # Set the coupling on channel 2 to DC
scope.write(":CHAN3:SCAL 0.2") # Set the scale to 0.2 V per division (square)

# Set the specifics for triggering a measurement
scope.write(":TRIG:EDG:SOUR CHAN2") # The measurement will use the signal coming in on channel 2 to trigger the scope
scope.write(":TRIG:EDG:LEV 0.1") # A measurement is triggered when the signal hits 0.1 V
scope.write(":TRIG:SWE NORM") # set the sweep to normal so the image is stable

# Set the horizontal scale for the measurement
""" 
    Setting the timebase on the scope to 20 ms gives a RAW sample rate of 125 MHz. 
    This allows for the entire signal to be seen on the screen of the scope. Only 
    1200 of the 3 million data points are displayed on the screen. When downloaded
    in the RAW state you will receive a large number of points in between the ones 
    shown on the screen allowing for much higher frequency resolution in the 
    spectrum. Using this setting allows easy acces to spectra with 10 Hz resolution.
"""
scope.write(":TIM:MAIN:SCAL 0.02") # the timebase will use 20 ms/division

#set the number of signals to average in the measurement
#scope.write(":ACQ:TYPE AVER")
scope.write(":ACQ:AVER 1")
scope.write(":RUN")

6


"""
    Use f_sour as an initial source frequency. This is used to verify the communication with
    fucntion generator.
    
    The function set_sour_f will set both channels on the function generator to supply a 
    900 Hz sine wave and then turn both outputs on. It will also make sure that the scope
    is in RUN mode and collecting data.
"""
f_sour = 900
#init = True # This is the initial setup of the function generator
do_spec_2.set_sour_f(scope,funcgen,f_sour, True)
time.sleep(5)


## This cell will read in the scope setup for data collection

"""scope.write(":STOP") # stop the scope
scope.write(":WAV:SOUR CHAN1") # set channel 1 as the source for data
scope.write(":WAV:MODE RAW") # use the raw data type to pull from internal memory
scope.write(":WAV:FORM ASC") # access the data in ascii format

# Read the preamble for the data download, to set physical parameters.
scope.write(":WAV:PRE?")
pream_dat = scope.read()
pream_dat = pream_dat.strip()
pream_dat = pream_dat.split(',')
float_pream_dat = [float(item) for item in pream_dat]
"""

float_pream_dat = do_spec_2.getpream(scope)

# Determine the Nyquist freqency for this measurement.
f_samp = 1/float_pream_dat[4] # Calculate the sampling frequency from dt
f_ny = f_samp/2 # The nyquist f is 1/2 teh sampling f
spec_res = 5.0 # The desired frequency resolution in Hz
n = int(f_samp/spec_res) # Figure the number of points needed to get a particular f resolution

print(f_samp, f_ny, spec_res, n)
init = False # No need to do the initial setup of FG

tb = 0.02 # timebase on the scope

for f in range(int(fmin), int(fmax) + int(df), int(df)):
    # Start a loop to rerun the iteration if fout is zero
    while True:
        scope.write(":ACQ:TYPE NORM")
        do_spec_2.set_sour_f(scope, funcgen, f, init)
        vscale_scope(scope, "CHAN1", tb)
        vscale_scope(scope, "CHAN3", tb)
        scope.write(":RUN")
        print("New frequency is ", f, " Hz")
        
        fout = do_spec_2.f_spec(scope, [1, 3], float_pream_dat[4], spec_res, f, n)
        
        # Check if fout contains zero values
        if (abs(value) == 0 for value in [fout[0].real, fout[0].imag, fout[1].real, fout[1].imag]):
            # If fout has valid (non-zero) values, break the loop and continue with the next part
            break
        else:
            # If fout has zero values, print a message and rerun the iteration
            print("Zero values detected, rerunning this iteration...")
    
    # If valid values are obtained, proceed
    print(fout)
    print(fout[0].real, "\t ", fout[0].imag, "\t", fout[1].real, "\t", fout[1].imag, "\n")
    fdat.write("%6.0f \t %5.2f \t %5.2f \t %5.2f \t %5.2f\n" % (f, fout[0].real, fout[0].imag, fout[1].real, fout[1].imag))

'''for f in range(int(fmin), int(fmax)+int(df), int(df)):
    scope.write(":ACQ:TYPE NORM")
    do_spec_2.set_sour_f(scope,funcgen,f,init)
    vscale_scope(scope,"CHAN1",tb)
    vscale_scope(scope,"CHAN3",tb)
    #do_spec_2.setVscale(scope,funcgen)
    #scope.write(":ACQ:TYPE AVER")
    #scope.write(ACQ:TYPE:AVER")
    scope.write(":RUN")
    #time.sleep(18)
    print("New frequency is ", f, " Hz")

    fout = do_spec_2.f_spec(scope,[1,3],float_pream_dat[4],spec_res,f,n)
    print(fout)
    print(fout[0].real, "\t ", fout[0].imag, "\t", fout[1].real ,"\t", fout[1].imag, "\n")
    fdat.write("%6.0f \t %5.2f \t %5.2f \t %5.2f \t %5.2f\n" %(f, fout[0].real,fout[0].imag, fout[1].real,fout[1].imag))'''

12500000.0 6250000.0 5.0 2500000
New frequency is  100  Hz

The index for the measurement frequency:  20
This is 100.000000 Hz

Magnitude at measurement frequency:
Chan 1: 7.212e+04
Chan 3: 7.663e+04

The maximum values and their indexes:
Chan 1: 7.212e+04, 20
Chan 3: 7.663e+04, 20
Chan1 data points;  2500000
Chan3 data points;  2500000
[(-66295.91879551366-28389.010329715973j), (-72107.51797987752-25947.645200738603j)]
-66295.91879551366 	  -28389.010329715973 	 -72107.51797987752 	 -25947.645200738603 

New frequency is  110  Hz

The index for the measurement frequency:  22
This is 110.000000 Hz

Magnitude at measurement frequency:
Chan 1: 7.935e+04
Chan 3: 8.368e+04

The maximum values and their indexes:
Chan 1: 7.935e+04, 22
Chan 3: 8.368e+04, 22
Chan1 data points;  2500000
Chan3 data points;  2500000
[(-58874.18420121051+53198.20162864868j), (-58523.89127314754+59816.50527865792j)]
-58874.18420121051 	  53198.20162864868 	 -58523.89127314754 	 59816.50527865792 

New frequency is  120  Hz

The index for the measurement frequency:  24
This is 120.000000 Hz

Magnitude at measurement frequency:
Chan 1: 8.273e+04
Chan 3: 8.766e+04

The maximum values and their indexes:
Chan 1: 8.273e+04, 24
Chan 3: 8.766e+04, 24
Chan1 data points;  2500000
Chan3 data points;  2500000
[(22959.53111038956+79476.21141743424j), (29624.90394263319+82504.78233084198j)]
22959.53111038956 	  79476.21141743424 	 29624.90394263319 	 82504.78233084198 

New frequency is  130  Hz

The index for the measurement frequency:  26
This is 130.000000 Hz

Magnitude at measurement frequency:
Chan 1: 8.569e+04
Chan 3: 9.083e+04

The maximum values and their indexes:
Chan 1: 8.569e+04, 26
Chan 3: 9.083e+04, 26
Chan1 data points;  2500000
Chan3 data points;  2500000
[(83313.54406396413+20026.329202841713j), (89537.89031629139+15266.430687252752j)]
83313.54406396413 	  20026.329202841713 	 89537.89031629139 	 15266.430687252752 

New frequency is  140  Hz

The index for the measurement frequency:  28
This is 140.000000 Hz

Magnitude at measurement frequency:
Chan 1: 9.233e+04
Chan 3: 9.827e+04

The maximum values and their indexes:
Chan 1: 9.233e+04, 28
Chan 3: 9.827e+04, 28
Chan1 data points;  2500000
Chan3 data points;  2500000
[(89817.00583023166+21399.505076229903j), (96845.13168084621+16686.0940606131j)]
89817.00583023166 	  21399.505076229903 	 96845.13168084621 	 16686.0940606131 

New frequency is  150  Hz

The index for the measurement frequency:  30
This is 150.000000 Hz

Magnitude at measurement frequency:
Chan 1: 1.034e+05
Chan 3: 1.110e+05

The maximum values and their indexes:
Chan 1: 1.034e+05, 30
Chan 3: 1.110e+05, 30
Chan1 data points;  2500000
Chan3 data points;  2500000
[(77434.79717186047-68530.10156721213j), (77997.4243739745-79009.88941879648j)]
77434.79717186047 	  -68530.10156721213 	 77997.4243739745 	 -79009.88941879648 

New frequency is  160  Hz

The index for the measurement frequency:  32
This is 160.000000 Hz

Magnitude at measurement frequency:
Chan 1: 1.202e+05
Chan 3: 1.298e+05

The maximum values and their indexes:
Chan 1: 1.202e+05, 32
Chan 3: 1.298e+05, 32
Chan1 data points;  2500000
Chan3 data points;  2500000
[(-17371.231954160532-118911.9765391805j), (-27380.76744584486-126836.59405257179j)]
-17371.231954160532 	  -118911.9765391805 	 -27380.76744584486 	 -126836.59405257179 

New frequency is  170  Hz

The index for the measurement frequency:  34
This is 170.000000 Hz

Magnitude at measurement frequency:
Chan 1: 1.481e+05
Chan 3: 1.612e+05

The maximum values and their indexes:
Chan 1: 1.481e+05, 34
Chan 3: 1.612e+05, 34
Chan1 data points;  2500000
Chan3 data points;  2500000
[(-136233.29256669706-57997.48137584438j), (-152207.89890476677-53157.53055499521j)]
-136233.29256669706 	  -57997.48137584438 	 -152207.89890476677 	 -53157.53055499521 

New frequency is  180  Hz

The index for the measurement frequency:  36
This is 180.000000 Hz

Magnitude at measurement frequency:
Chan 1: 1.523e+05
Chan 3: 1.634e+05

The maximum values and their indexes:
Chan 1: 1.523e+05, 36
Chan 3: 1.634e+05, 36
Chan1 data points;  2500000
Chan3 data points;  2500000
[(-127135.49524085893+83807.96754080558j), (-130270.48146198326+98687.9699401313j)]
-127135.49524085893 	  83807.96754080558 	 -130270.48146198326 	 98687.9699401313 

New frequency is  190  Hz


# Run this cell to make sure to close down the connection with the USB ports correctly. 
# If this cell is not run. Communications will be disrupted.
funcgen.write(":OUTP1 OFF") # Turn off the output of channel 1 (Don't forget this step or nothing will happen)
funcgen.write(":OUTP2 OFF") # Turn on the output of channel 1 (Don't forget this step or nothing will happen)

scope.write(":CHAN1:DISP OFF") # Turn off the display of channel 1
scope.write(":CHAN2:DISP OFF") # Turn off the display of channel 2
scope.write(":CHAN3:DISP OFF") # Turn off the display of channel 3

funcgen.close()
scope.close()
fdat.close()
print('Ended program ready to run again')

Ended program ready to run again


#write ab if mic A at 11.5, B at 1.5, or ba if ic A at 1.5, B at 11.5
"""
fdat_ab = open(path_to_file, "r")
fname = 'T_ab_'+str(td.year)+"-"+str(td.month)+"-"+str(td.day)+"-"+ts.strftime("%H%M%S")
path_to_tfile = fpath + fname + ".dat"
file_exists = exists(path_to_tfile)
if file_exists:
    print("old file exists, rename transfer output file")
else:
    tranx_file = open(path_to_tfile, "w")
TFC.Transfer_Initial(fdat_ab, tranx_file)
tranx_file.close()
"""

This note book will be used to record impedance tube data for the A50-150 plug.¶